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ABSTRACT 

We investigate the question of the proper thermal averaging of neutralino annihilation 
amplitudes which possess poles and thresholds, as they impact on the calculated neutralino 
relic density and therefore on the cosmological viability of supersymmetric models. We 
focus on two typical resonances, namely the Z boson and the lightest Higgs boson (h). 
In the context of supergravity models with radiative electroweak symmetry breaking, an 
exploration of the whole parameter space of the model is possible and the overall relevance 
of these sophisticated analyses can be ascertained. As an example we chose the minimal 
SU (5) supergravity model since the presence of such poles is essential to obtain a cosmo- 
logically acceptable model. We find that the proper thermal averaging is important for 
individual points in parameter space and that the fraction of cosmologically acceptable 
points is increased somewhat by the accurate procedure. However, qualitatively the new 
set of points is very similar to that obtained previously using the usual series approxima- 
tions to the thermal average. We conclude that all phenomenological analyses based on 
the previously determined cosmologically allowed set remain valid. 
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1. Introduction 



Much speculation has gone on for some time about the nature of the "observed" 
astrophysical dark matter in the Universe |jj . The lightest supersymmetric particle - the 
lightest neutralino - is a prime candidate for a cold dark matter relic @|§, and as such 
it would constitute an essential ingredient in contemporary structure formation ideas f|]. 
A consistency check for a possible dark matter candidate x i n the Big-Bang cosmology is 
provided by an independent lower bound on the age of the Universe. The present relic 
abundance fi x = p x / p cr it must be bounded above by n x /ig < 1 f|, where 0.5 < ho < 1 
is the Hubble parameter in units of 100 kms -1 Mpc -1 . Weaker bounds follow from direct 
astrophysical determinations of Oq- Also, in inflationary cosmology Oq = 1 an d the 
observational bound is also satisfied. 

In the specific case of the lightest neutralino, the computation of Q x has been at- 
tempted in several supersymmetric (minimal 0,^[7|, non-minimal H, supergravity with- 
out PJim, and with [|II| , [r2|,|IB| radiative electroweak breaking) models to varying degrees 
of approximation for nearly a decade. The details of the particle physics model come into 
play mainly in the calculation of the total annihilation amplitude XX ~ y a ^ an d its thermal 
average. Under normal circumstances, thermal considerations allow one to conclude that 
center-of-mass energies close to the minimal one (-^/i = 2m x ) are the most likely ones H, 
and that a rapidly converging series expansion around this point should suffice ||14|| . For 
two-body final states with masses ra\^ to be "open" (in the non-relativistic limit) it is 
necessary that m\ + rri2 < 2m x . Since x is the lightest supersymmetric particle, only non- 
supersymmetric final states contribute (i.e., qq, Z + Z~, W + W~ , ZZ, hh, hA, hH, Zh, ZA, ■ ■ -, 
where h, A, H are the supersymmetric Higgs bosons). 

Such series expansions, however, have been noted to fail badly when s-channel reso- 
nances and/or new-channel thresholds are present in the annihilation amplitude |T5,rB[. 



The discussion in these references provided a solid ground for tackling such problems, but 
their practical applications were only briefly explored. On the other hand, supergravity 
models with radiative electroweak symmetry breaking |17| provide a fertile testing ground 
for such sophisticated techniques, since the complete mass spectrum and couplings of the 
model can be specified in terms of very few parameters |JL8|| . In particular, the occurence of 
poles and thresholds in the annihilation amplitude can be studied over the whole parameter 
space of these models, in order to determine whether the more accurate results provided by 
these techniques change significantly the cosmologically allowed region of parameter space 
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or not. These special considerations play a very important role in the cosmo logical via- 
bility of the longest-standing GUT model, namely the minimal SU(5) supergravity model 
fl9|] . It has been recently pointed out |1^,|13] that, because the stringent proton decay con- 
straints on this model force the most efficient annihilation channel mediators to be very 
heavy, the neutralino relic density is small enough basically only near the lightest Higgs 
(h) and Z-boson resonances, i.e., for m x ps hrrih^z- Thus, a more accurate treatment, as 



described in Refs. |T^,^6[, appears mandatory in this case. 

The purpose of this paper is twofold. Firstly, we present a quantitative discussion 
of the accurate thermal averaging necessary in the context of this class of models, and 
explicitly show the breakdown of the usual series expansion near the poles and its not- 
so-good accuracy not-so-near the poles. Secondly, we perform a complete recalculation 
of the relic density in the minimal SU(5) supergravity model, and show the effects of 
the more accurate treatment and the explicit role of the s-channel resonances (h, Z) and 
the XX ~~ > hh threshold. We find that even though the re-computed values of the relic 
density are shifted relative to our previous (less accurate) results, the overall fraction of 
the parameter space which is cosmologically allowed is not qualitatively changed, and thus 
all predictions based on the previously determined cosmologically allowed set remain valid. 



2. The thermal average 

The neutralino relic density is given by the expression WA 



VL x h 2 = 1.555 x 10 8 (m x / GeV)fc(0)g(0), (2.1) 

where h(0) = 3.91 is the effective number of entropy degrees of freedom today, and q(0) 
is obtained by solving the Boltzmann equation for q = n/(T 3 h(T)), with n the actual 
number density of x particles at temperature T. This equation is given by 

^ = A(x)( ? 2 - q%)(x) (2.2) 

with x = T/m x , qo the analog of q but with the x particles in thermal equilibrium, and 

X(x) = {±n 3 G N y 1/2 [h(T) + \m x xh'{T)] (av m ), (2.3) 

where Gjy is the gravitational constant and g(T) the effective number of energy density 
degrees of freedom. For a detailed discussion of how to evaluate all the terms appearing 
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in the Boltzmann equation and how to solve the equation itself, see e.g., the Appendix 



in Ref. [^0|. The novelty in the present discussion is in the evaluation of the thermal 
average factor (crvMfii), where vm£\ is the M0ller velocity |T5[]. The general expression for 
this quantity is [[15] 

1 f°° 

{aVM * l) = Sm+TKUx-^ L> a ' (S ~ 4 -x)v^^i(v^/T) da, (2.4) 

where Ki are the modified Bessel functions of order i, and cr is the total annihilation cross 
section. To make contact with previous work we rewrite {(JVm^i) in Lorentz invariant form 
as follows 

1 f°° 

((niMfii) = 4 m 5 xK 2( x -i^ J 4 2 ds ( s - ^rn 2 x ) 1/2 K x {yfs/xm x ) w(s), (2.5) 



with [H| 

w(s) = \ f dLIPS \A{xX - *U)\ 2 - (2-6) 

The usual series expansion for the thermal average follows from the observation that only 
for x ^ 0.1 is the Boltzmann equation sensitive to the value of A(x). In this regime the 
argument of K\ in Eq. (|2.5[) is always larger than ^/s/xm x > 2/x ^ 20, and the Bessel 



function (Ki(y) ~ ^n/2ye~ v , y ^> 1) dies away quickly with increasing ^fs. Therefore a 
series expansion of w(s) around i/i = 2m x should converge quickly. The resulting series 
of integrals can be done analytically giving [|T4| 



(cn^i) = —[«?- §(2w - w')x + §(16«; - 8w' + 5w")x 2 + C(a; 3 )] s= 4 m 2 (2.7a) 
m x 

= a + bx + cx 2 + 0(x s ). (2.7b) 

The problem with this approximation when w(s) has a pole can be best seen in a 
simple analytical example. Let us consider the case where an s-channel resonance with 
mass m R and width F R dominates w(s), i.e., 

y(y - 1) 
(y - vr) 2 + 1r 

where y = s/Am 2 x , y R = m 2 R /Am 2 x , and ^ R = Y R m R /Am 2 x = (V R /m R )y R . This form in 
fact applies to XX Z — > / / when m/ = 0, with C some dimensionless function of the 
couplings (see Eqs. A. 2 and A. 3). We then obtain 

2C f» y(y-l) 



w(y)=C - y \\, [ 2 , (2.J 



2/(2/-!) /l^l.-a^-D/a 



<fo , ^ s « ^ e -*Wv-V/*, (2.9) 



^rm^x 3 / 2 A (?/ - y R ) 2 + 7I V v 7 ^ 



where the limiting form holds in the x-regime of interest (x ^ 0.1). Expanding w(y) around 
y = 1 and performing the integrals analytically one arrives at the series expansion in Eq. 
(fT7p with 



0, (2.10a) 
It- ~Trt 2-, (2.106) 

f l - . 10(1 7 ^ } . 1 • (2.10c) 
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(i-y*) 2 + 7| 



For illustrative purposes we have computed (crv M ^{} exactly using Eq. (2.9) and also using 
the expansions in Eqs. (|2.7|) ) (|2.10|) . The results are shown in Fig. 1 (2) (in arbitrary 
units) for the Z (/i)-pole, where we have taken x = 0.05 and Tr/uir = 0.0274(0.0002) 
for R = Z (h). Clearly, sufficiently away from the poles (yR > 1 or t/fl < 1) the series 
expansions approach the exact result, while near the poles these are highly innacurate. In 
particular, the second-order approximation (a + bx + ex 2 ) fails badly right above the poles 
(yR < 1) since it gives negative thermal averaged annihilation cross sections. 

The degree of inaccuracy of the series espansions can be better appreciated by exam- 
ining the ratios (cn?M«si) exact/ (ct'M^i)approx, as shown in Fig. 3 for both resonances. Below 
the pole (yR > 1) one can clearly see the quicker convergence of the higher-order approxi- 
mation. The peak in the a + bx + ex 2 line is due to the series expansion changing sign at 
y = 0.75, giving negative (non-sensical) results for 0.75 < y & 1. The overall result is that 
as one approaches the poles from below (yR > 1) the exact thermal average first becomes 
larger than naively expected (relic density smaller) , then near the poles it becomes smaller 
(relic density larger), until above the poles where it quickly approaches the naive estimate 
(in first order, a + bx). In practice this means that the neutralino relic density distributions 
as a function of the neutralino mass will show not-as-narrow (broader) and not-as-deep 
(shallower) pole structures, and that these will be asymmetric. 

The bottom row in Fig. 3 shows a detail of the convergence of the series expansions 
to the exact result below the pole. It is perhaps somewhat unexpected that the expansions 
are relatively innacurate not-so-near the poles, i.e., for yR = 2, m x « 0.7(^m^), which 
for the Z-pole gives m x « 32 GeV. The reason for this behavior can be understood 
by studying the integrand in Eq. (2.9) with w(y) in exact and Taylor-expanded (around 
y = 1) forms. These are shown as functions of y in Fig. 4 for both poles and representative 
values of yR. The solid/dotted/dashed lines correspond to the integrand evaluated using 
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w(y) in exact /up-to-first-derivates/up-to-second-derivatives forms. Clearly, the poles at 
V = Ur are 'missed' by the expansions for y^ <; 1. However, for sufficiently high values of 
yu the exponential kills-off the poles completely, although they 'last longer' for the much 
narrower /i-pole. Note also the quantitative effect of the expansion up to second-derivatives 
(dashed lines) relative to that up to first-derivatives (dotted lines). In practical relic 
density calculations only the first derivatives of w(y) are usually kept (which corresponds 
approximately to (ctvm^i) = a + bx) and therefore good accuracy is not reached until 
considerably away from the poles (from below; from above convergence is faster). Since 
w(y) usually receives additional contributions from non-resonant channels which are likely 
to overshadow the contributions from resonant channels when away from the corresponding 
resonances, in practive the behavior away from the poles is not easily distinguishable. 



3. Application to the minimal SU(5) supergravity model 

As mentioned above, the study of the relic density of neutralinos requires the knowl- 
edge of the total annihilation amplitude XX ~^ The latter depends on the model 
parameters to determine all masses and couplings. Previously [!(]] we have advocated the 
study of this problem in the context of supergravity models with radiative electroweak 
symmetry breaking, since then only a few parameters (five or less) are needed to specify 
the model completely. In particular, one can explore the whole parameter space and draw 
conclusions about a complete class of models. The ensuing relationships among the various 
masses and couplings have been found to yield results which depart from the conventional 
minimal supersymmetric standard model (MSSM) lore, where no such relations exist. Here 
we study a novel aspect of these correlations, namely the ocurrence of poles and thresholds 
in the annihilation amplitude. We choose to work with the minimal SU(5) supergravity 



model |nj since its five-dimensional parameter space is strongly constrained by the proton 
lifetime and a sufficiently small neutralino relic density. In fact, the latter is cosmolog- 
ically acceptable only because of enhancements in the neutralino annihilation amplitude 
near s-channel Z and h resonances. Our purpose here is to determine whether an accu- 
rate computation of the thermal avarage changes significantly the previously (inaccurately) 
determined cosmologically allowed region of parameter space. 

We have performed an extensive search of the five-dimensional parameter space of 
the model. The five parameters are: the top-quark mass (mt), the ratio of Higgs vacuum 
expectation values (tan/3), and three universal soft-supersymmetry breaking terms (the 
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gaugino mass m^, the scalar mass mo, and the trilinear scalar coupling A). The sign 
of the Higgs mixing parameter \x is also undetermined. Our search accepts only those 
points which give adequate radiative electroweak symmetry breaking and satisfy all known 



phenomenological constraints in the sparticle and Higgs spectrum, as described in Ref. JT^] . 
We also include the very restrictive proton decay constraint pU|JT3|] with the unification 



scale calculated using two-loop gauge coupling unification including the effect of light 
supersymmetric thresholds ||21|| . The remaining points in parameter space (~ 2000 per 
sign of fx), which have tan/3 = 1.5,1.75,2.0, were then used to compute the neutralino 



relic density using the methods of Ref. |10| for solving the Boltzmann equation and the 
accurate thermal averaging procedure of Ref. |T3[, as described in Sec. 2. 

Since in this model m x ^ 65 GeV (and even m x ^ 50 GeV if improved lower bounds on 
rrih are further imposed [p2|), and the Higgs masses obey rrih ^ 100 GeV, ttia,h > 500 GeV, 
only the h and Z s-channel resonances contribute to XX ~ * ffi with / a light quark or 
lepton.0 In addition, the channel XX ~ > hh can be kinematically allowed (showing a 
threshold effect) in some regions of parameter space (m x > rrih). Note that in the present 
approach to calculating (ctvm^i), in principle all kinematically allowed channels contribute, 
since y/s > 2m x is allowed. This is in contrast with the traditional approach where 
y/s = 2m x is fixed and only two-body channels with rrii + rri2 < 2m x contribute. However, 
the farther away y/s is from its lower limit of integration in Eq. (|2.5|) , the least it will 
contribute to the total integral (see Eq. (2.9) where y/y = y/s/2m x ). This is exemplified in 
Fig. 4 which shows that for y <; 2, the contribution to the integral is negligible. Important 
annihilation channels that one would need to worry about (besides XX ffi hh) include 
XX ~^ W + W~ , ZZ, which have /i-pole and Z-pole (only VF + VF~) contributions when the 
neutralino is not a pure gauge eigenstate. For m x ^ 50 (65) GeV, these channels open up 
for y 2.6,3.4(1.55,2.0) respectively. Thus, only the W + W~ channel (which requires 
y/s > 2Mw) could be relevant. To give any significant contributions, i.e., y = s/4m x ^ 2, 
one would need m x ^ 57 GeV, but in this case yz < 0.65 and yn < (m^/114 GeV) 2 ^ 0.77 
(for rrih & 100 GeV), and the poles are not encountered. Therefore, in the following 
we have only included the XX ~^ //> hh annihilation channels. The exact expressions 
for w(s) are given in the Appendices, including all interference terms, which have been 
typically neglected in previous analyses. We note that even near the poles (for y ^ 3) the 



1 In our calculations, the /i-width has been obtained for every point in parameter space includ- 
ing the h — > ff and h — > X1X1 (when kinematically allowed) channels. 
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interference terms are generally smaller than the squared terms, although for large y (when 
all contributions are unimportant anyway) the two contributions can be comparable. 

The results of our calculations are shown in Fig. 5 and Fig. 6 where (ctvm^i) has been 
computed exactly and in the series expansion to first order (a + bx), respectively. The 
overall result is not very clear from these figures, although some shifts of the points are 
evident. In particular, around the Z-pole (m x ~ \Mz) one can see that the pole structure 
in the exact case is broader, shallower, and asymmetric (lower below the pole) relative 
to the approximate solution, as anticipated in Sec. 2. Moreover, there is a ~ 53% (27%) 
increase in the number of cosmologically allowed points (fi x /io < 1) for // > (/i < 0) 
relative to the approximate result. This shift is however not qualitatively significant. In 
fact, the distribution of cosmologically allowed points in the (m ±,m^) plane (see Fig. 
7) is very close to that obtained previously using the series expansion (see Fig. 2 in Ref. 

To show the effect on the relic density of poles and thresholds of the annihilation 
amplitude, we show in Fig. 8 in the (m x , rrih) plane those points in parameter space where 
Qexact i^approx > 1 These points correspond to ((rv M ^) exact / (av M ^) approx < 1, and 
according to Sec. 2 (see Fig. 3) should occur very near and above the poles (yz,h !)• This 
is precisely what Fig. 8 shows, where the Z-pole (m x ~ \Mz) and the h-pole (m x ~ \rrih) 
are denoted by dashed lines. In Fig. 9 we show those points where Q^ xact /Qapprox <- q g 
which correspond to (vvmai) exact / (gvmm) approx ^ 1) and according to Fig. 3 should occur 
below the poles (yz,h *5 !)■ This is again borne out by the results in Fig. 9. Moreover, 
for m x ~ rrih the XX ~^ hh channel opens up and a threshold effect causes a drop in 
(O x /iQ) exact . These points show up in Fig. 9 along the diagonal. It is interesting to 
note that these points lie sligthly above the diagonal, i.e., m x ^ rrih. In the usual series 
approximation, this channel would not open up until m x > rrih. I n the exact treatment, 
even for m x ^ rrih, \fs > 2rrih is possible and the channel becomes kinematically allowed. 
However, as discussed above, only values of \fs close to its minimum value (y/s = 2m x ) 
can contribute significantly, and this is why the effect only occurs very near the diagonal. 
For m x > rrih both methods give similar results. 

4. Conclusions 

We have investigated the question of the proper thermal averaging of neutralino an- 
nihilation amplitudes which possess poles and thresholds, following the methods of Ref. 
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T5f| . We have focused on two typical resonances in supersymmetric models, namely the 



Z boson and the lightest Higgs boson (h). In the context of supergravity models with 
radiative electroweak symmetry breaking, an exploration of the whole parameter space 
of the model is possible and the overall relevance of these sophisticated analyses can be 
ascertained. As an example we chose the minimal SU(5) supergravity model since the 
presence of such poles is essential to obtain a cosmologically acceptable model. We have 
found that the proper thermal averaging is important for individual points in parameter 
space. Also, the fraction of cosmologically acceptable points is increased somewhat by 
the accurate procedure. However, qualitatively the new set of points is very similar to 
the previously allowed set.! We conclude that all phenomenological analyses based on the 
previously determined cosmologically allowed set remain valid. 
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Concurrently with our calculation there has appeared an analogous one [23| which reaches 
similar qualitative conclusions. More quantitative comparisons are not possible given the lack of 
detail in Ref. ||§. 
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Appendix A. Function w(s) 

In this Appendix we present the full set of explicit expressions for the Lorentz invariant 
function w(s) that we used in this work. In what follows, all the couplings except those 
specified otherwise can be found in Ref. [[nj. As discussed in Sec. 3, in the minimal SU (5) 
supergravity model, the only relevant annihilation channels are XX ~* ff an d XX ~^ hh, 
therefore, w(s) can be written as 



w(s) 



32tt 



Am 



^c f B{s - 4m}) m f -w {ff \s) + 9{s - 4mjWl ±w {hh \s) 



Ami 



(A.l) 

where Cf is the color factor of conventional fermion / (c/ = 3 for quarks and Cf = 1 for 
leptons), the summation in the first term runs over all conventional quarks and leptons 
except the top quark. 

A.l. XX ~ * // annihilation 

There are totally eight types of contributions to w^^(s): 
(1) Z-boson exchange: 



\G 



xxZ\2 ( m 2 m 2 

} J ioir //Z l 2 x f (» M 2 \ 2 

-3( s _Ml)2 + Tpf|1^ |GA 1 ~M% KS Mz) 



+ 



\G f A fz \ 2 (s - Am}) + \G f v fz \ 2 (s + 2m})] (s - Am 



i) 



(A.2) 



(2) CP-even Higgs-boson (Si = h, S2 = H) exchange: 



w { / f \s) 



r<xxSi^if fSi 



(3) CP-odd Higgs-boson (A) exchange: 



(s - Am})(s - 4m£); 



(A.3) 



w { l f \s) 



r<xxA n ffA 

p p 



s — m 2 A + iT aw-a 



(A.4) 



(4) Z-A interference: 



w 



(If) 

ZP 



00 



8em x mf 



M 2 



Re 



r<xxA n ffA 

p p 



^A ^A 



s-m A + iF A m A I \ s - Mf + iT z M z 



(M 2 z -s)s; (A.5) 
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(5) sfermion (/) exchange: 



w 



iff) 



1,3 



T 2 - 2(m 2 x + m})T x + ({m 2 x + m)f + 4m 2 x m 2 f )T 

+ (AiAi - B l _B J _) |T 2 - 2(m\ + m))T x + {m\ - m 2 f ) 2 T 
+A\A J + (m 2 x + m 2 )s-Am 2 x m 2 y 1 + A l _A°_ {m\ - mj)s^i 
+ (B l + B + ~ B - B -) [^2 + (ml + m)) 2 y x 
+ (A%B J + + A\B\)em x m s [4(m^ + m 2 f )T - AT X + s^i] 

-{A\B j + - A{B l + )em x m f y Ms, mi 2 ,m 2 f ,m 2 ~ ,m 2 ~); (A.6) 

where T 2 , %, y 2 , 3^i, 3^o are some auxiliary functions which are given in Appendix B, 
and the summation runs over two sfermion mass eigenstates fi j2 for each corresponding 
fermion /. Also, we have introduced the following notation for sfermion couplings: 



B\ 



with 



G i =G* ffL U u + G 



xffn 



u 2i , 



G R =G R ffL U ll + G R JJ «U 2 



xffi 



ii 



(A.7a) 
(A.76) 

(A.8a) 
(A.86) 



where U is the orthogonal matrix which rotates /l,r into fi t2 . 
(6) Z-f interference: 



w {I l\s) = Re 



r<xxz s-iffZ 

A A \ 

s-Ml + iT z M z / zf 



+ Re 



^ a b y 
s -Ml + iT z M z 



I {2 1 
zf 



(A.9) 



where 

^7 = 5> S WW ( M2 



3 )~ ( s + (l72 - 3 )( m x + m 2 f -m 2 j))F(s,m 2 x ,m 2 f ,m 2 j) 



M 2 



8m 2 m 2 f 



-A\ \(2{m 2 + m 2 f )s H tpt^-s ~ 2(m£ + m 2 f - m~) 2 - 16m 2 m 2 f ) T{s, m 2 m 2 f , m~) 



+ s + 2(m 2 x + m 2 f -m 2 -) 



(A.lOa) 



J z~ = S A ~ ( 2 ( m i ~ m )) s - 2 ( m x + m 2 f -m 2 y) 2 + 8m 2 x m)) F(s, m 2 , m 2 f , m 2 j) 

i 

+ s + 2(m x + m 2 - m 2 j) . (A. 106) 
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(7) S-f interference: 



w 



iff) 
sf 



w = E Re — 



rixxSj n ffSj 



i=i,2 



s - m s j + iT s j m Sj J sr 



(A.ll) 



where 



i" s ^=^< — 2A!^em x m/ 2 + (s — 2(m x + m 2 — mj))j r (s,m x ,m 2 ,m~) 



s + + m 2 — m~)s — 8m x m 2 )j-'(s 1 m x , m^, m~) >. (A. 12) 



(8) A-f interference: 



where 



w ~{s) = Re ' 
Af v ' 



s -m\ + zr A m A ) 



(A.13) 



1 ' aJ = E I ^A l + em x mfsJ r {s, m 2 x , m 2 , m y) 



+B l + -s + (m x + mj - m~)sJ r (s, m x , m 2 , m~) 



(A.14) 



A. 2. xx — > /i/i annihilation 

There are totally three types of contributions to w^ hh \s): 
(1) CP-even Higgs-boson (Si = h, S2 = H) exchange: 



riXXSi/-ihhSi 

ys u 

s — 



s ~ m s l + iT s i m Si 



(s - 4m x ); 



(A.15) 



(2) neutralino exchange: 



w 



<W( S ) = - \G^ h \ 2 \G x s x > h \ 2 { T 2 r 



s + 2(m 2 x - m 2 h ) + 2tm x (e % m Xi + e 3 m X] ) 



- \^itj m Xi m xA s - 4m x) - ( m x - m h) ( m l - m l + 2em x (e,m Xi + e 3 m Xj )) 



To 



+ 



(e i e J m Xi m Xj + em x (e;m Xi + ejm Xj )) (s - 4m x ) - [m 2 x - ml)(3m x + ml) 



yi 



+y2 + em x (e,m Xi - ej-m^XVoj (s, m£, m 2 h , m 2 Xi , m x .); 



(A.16) 
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(3) S-x inteference: 



where 



w 



i=i,2 



X' 



(A.17) 



J Sx = ^|G^ Xi/ *| 2 j2em x 1 + (m 2 + m Xi - m£).F(s,TO£,TO£,m Xi ) 

- e i TO Xi ( s - 4m x ) ^(s, m x , m|, m Xi ) J; (A.18) 



Appendix B. Some auxiliary functions 

In this Appendix we collect various auxiliary functions that appear in Appendix A. 
First, we define 



1 1 



t±{s, x,y) = x + y - -s± - y/ (s - 4x)(s - Ay) 



then 



and 



F(s,x,y,z) 



y/(s-4x){s-4y) 



In 



t+(s,x,y) - z 



t-(s,x,y) - z 



(B.l) 



(B.2) 



T 2 (s,x,y,z u z 2 ) = 1 + 
T 1 (s,x,y,z 1 ,z 2 ) 
T (s 7 x,y,z 1 ,z 2 ) 



yi(s,x,y,Z!,Z2) 
3>o(s,a,y,zi,z 2 ) 



{z\ - z 2 ) 



z\T{s, x, y, zi) - z\T{s, x, y, z 2 ) 



1 



{z\ - z 2 ) 

1 



{z\ - z 2 ) 
y 2 {s,x,y,z 1 ,z 2 ) = 1 + 



ziF{s, x, y, zi) - z 2 T{s, x, y, z 2 ) 

.F(s,a;,y,zi) - T{s,x,y,z 2 ) 
1 



(s - 2x - 2y + z\ + z 2 ) 

+ z 2 (s + z 2 -2x- 2y)J 7 (s, x, y, z 2 ) 
1 



(B.3o) 
(B.36) 
(B.3c) 

zi(s + zi-2x- 2y)F(s, x, y, zx) 

(B.3d) 



(s - 2x - 2y + z\ + z 2 ) 
1 



.F(s,a;,y,zi) + J 7 (s,x,y,z 2 ) 



(s-2(x + y- zij) F(s, x, y, z x ) 



(s - 2x - 2y + zi + z 2 ) 

-(s-2(x + y- z 2 )) .F(s, x, y, z 2 ) 



(B.3e) 



(B.3/) 
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Figure Captions 



Fig. 1. The thermal average {ctvm^i) (in arbitrary units) as a function of y z = 
(Mz /2m x ) 2 for a Z-pole dominance situation. The exact and traditional series 
expansion results are shown. Note that the second-order approximation gives 
negative thermal averaged annihilation cross sections right above the pole. 

Fig. 2. Same as Fig. 1 but for an /i-pole dominance situation. 

Fig. 3. The ratios (av m^i) exact / {cvm^i) a PP rox for the two approximate series expansions 
(first-order: dashed; second-order: solid) and for the two sample poles considered. 
The peak in the a + bx + cx 2 line is due to the series expansion changing sign 
at y = 0.75, giving negative (non-sensical) results for 0.75 < y ^ 1. The bottom 
row plots show the slow convergence not-so-near the poles to the exact result. 

Fig. 4. The integrand in the thermal average (Eq. (2.9)) as a function of y for two 
representative values of y^ and the two poles considered. The w(y) function is 
used in exact form (solid lines) and Taylor-expanded form up to first (dotted 
lines) and second (dashed lines) derivatives. Note that the expansions miss the 
poles completely, although these are killed-off when sufficiently away from the 
poles. 

Fig. 5. The calculated values of the neutralino relic density in the minimal SU(5) su- 
pergravity model as a function of m x using the exact thermal average procedure. 
Note that the pole structure near the Z-pole is broader and shallower than that 
shown in Fig. 6 where the series approximation to the thermal average has been 
used. 

Fig. 6. Same as Fig. 5 but using the first order series approximation to the thermal 

average ((cn>M^i) = a + bx). 

Fig. 7. Those points in parameter space where ^ x h,Q < 1 with the thermal average com- 
puted exactly. These results do not differ qualitatively from those shown in Fig. 
2 in Ref. where the series approximation to the thermal avarage was used. 



Fig. 8. Those points in parameter space where Q^ xact /^approx > ^ which according to 
Fig. 3 should occur near and above the poles. The dashed lines indicate the Z- 
and /i-poles. 

Fig. 9. Those points in parameter space where Q^ act j 'Qapprox ^ q w ] 1 j c ] 1 according 
to Fig. 3 should occur below the poles. The dashed lines indicate the Z- and 
/i-poles, as well as the XX ~^ hh threshold. 
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